% 智能电网调度优化测试脚本
% 基于SAEA-RFS算法的电网调度系统测试

clear;
clc;
close all;

% 屏蔽RBF矩阵奇异值警告
warning('off', 'MATLAB:nearlySingularMatrix');
warning('off', 'MATLAB:singularMatrix');

fprintf('========================================\n');
fprintf('智能电网调度优化系统测试\n');
fprintf('基于SAEA-RFS算法\n');
fprintf('========================================\n\n');

%% 1. 初始化电网模型
fprintf('1. 初始化电网模型...\n');
grid_data = grid_model();

%% 2. 设置算法参数
fprintf('\n2. 设置算法参数...\n');
algorithm_params = struct();
algorithm_params.population_size = 15;      % 种群大小
algorithm_params.initial_samples = 80;      % 初始训练样本数量
algorithm_params.sub_problems = 12;         % 子问题数量
algorithm_params.sub_problem_size = 70;     % 子问题最大维度
algorithm_params.max_iterations = 8;        % 子问题优化最大迭代次数
algorithm_params.max_evaluations = 3000;    % 最大函数评估次数

fprintf('算法参数设置完成：\n');
fprintf('- 种群大小: %d\n', algorithm_params.population_size);
fprintf('- 初始样本数: %d\n', algorithm_params.initial_samples);
fprintf('- 子问题数量: %d\n', algorithm_params.sub_problems);
fprintf('- 子问题维度: %d\n', algorithm_params.sub_problem_size);
fprintf('- 最大迭代次数: %d\n', algorithm_params.max_iterations);
fprintf('- 最大评估次数: %d\n', algorithm_params.max_evaluations);

%% 3. 运行SAEA-RFS优化算法
fprintf('\n3. 开始电网调度优化...\n');
[best_schedule, optimization_history, total_time] = grid_saea_rfs(grid_data, algorithm_params);

%% 4. 结果分析
fprintf('\n4. 结果分析...\n');

% 计算成本分解
[final_cost, cost_breakdown] = grid_objective(best_schedule.generator_power(:)', grid_data);

fprintf('\n最终调度方案成本分析：\n');
fprintf('- 总运行成本: %.2f 元\n', final_cost);
fprintf('- 发电成本: %.2f 元 (%.1f%%)\n', cost_breakdown.generation_cost, ...
    100*cost_breakdown.generation_cost/final_cost);
fprintf('- 网络损耗: %.2f 元 (%.1f%%)\n', cost_breakdown.network_loss, ...
    100*cost_breakdown.network_loss/final_cost);
fprintf('- 碳排放成本: %.2f 元 (%.1f%%)\n', cost_breakdown.carbon_cost, ...
    100*cost_breakdown.carbon_cost/final_cost);
fprintf('- 储能成本: %.2f 元 (%.1f%%)\n', cost_breakdown.storage_cost, ...
    100*cost_breakdown.storage_cost/final_cost);
fprintf('- 约束违反惩罚: %.2f 元 (%.1f%%)\n', cost_breakdown.violation_penalty, ...
    100*cost_breakdown.violation_penalty/final_cost);

%% 5. 可视化结果
fprintf('\n5. 生成可视化结果...\n');

% 创建图形窗口
figure('Position', [100, 100, 1400, 900]);

% 子图1：优化收敛曲线
subplot(2, 3, 1);
semilogy(optimization_history(:, 1), optimization_history(:, 2), 'b-', 'LineWidth', 2);
title('SAEA-RFS算法收敛曲线', 'FontSize', 12, 'FontWeight', 'bold');
xlabel('函数评估次数', 'FontSize', 10);
ylabel('总运行成本 (log scale)', 'FontSize', 10);
grid on;
legend(['最终成本: ', num2str(final_cost, '%.2f'), ' 元'], 'Location', 'best');

% 子图2：发电机出力曲线
subplot(2, 3, 2);
gen_power_plot = best_schedule.generator_power;
plot(1:24, gen_power_plot(1:min(5, size(gen_power_plot, 1)), :), 'LineWidth', 2);
title('发电机出力曲线 (前5台)', 'FontSize', 12, 'FontWeight', 'bold');
xlabel('时间 (小时)', 'FontSize', 10);
ylabel('出力 (MW)', 'FontSize', 10);
grid on;
legend(arrayfun(@(i) sprintf('发电机%d', i), 1:min(5, size(gen_power_plot, 1)), 'UniformOutput', false), 'Location', 'best');

% 子图3：储能设备充放电曲线
subplot(2, 3, 3);
storage_power_plot = best_schedule.storage_power;
plot(1:24, storage_power_plot(1:min(3, size(storage_power_plot, 1)), :), 'LineWidth', 2);
title('储能设备充放电曲线 (前3台)', 'FontSize', 12, 'FontWeight', 'bold');
xlabel('时间 (小时)', 'FontSize', 10);
ylabel('功率 (MW, 正值充电)', 'FontSize', 10);
grid on;
legend(arrayfun(@(i) sprintf('储能%d', i), 1:min(3, size(storage_power_plot, 1)), 'UniformOutput', false), 'Location', 'best');

% 子图4：负荷分配曲线
subplot(2, 3, 4);
load_allocation_plot = best_schedule.load_allocation;
plot(1:24, load_allocation_plot(1:min(5, size(load_allocation_plot, 1)), :), 'LineWidth', 2);
title('负荷分配系数 (前5个负荷)', 'FontSize', 12, 'FontWeight', 'bold');
xlabel('时间 (小时)', 'FontSize', 10);
ylabel('分配系数', 'FontSize', 10);
grid on;
legend(arrayfun(@(i) sprintf('负荷%d', i), 1:min(5, size(load_allocation_plot, 1)), 'UniformOutput', false), 'Location', 'best');

% 子图5：总发电量与负荷对比
subplot(2, 3, 5);
total_gen = sum(best_schedule.generator_power, 1);
total_load = sum(best_schedule.load_allocation .* ...
    arrayfun(@(i) grid_data.loads(i).base_load * grid_data.loads(i).load_profile, 1:grid_data.num_loads, 'UniformOutput', false), 1);
total_storage = sum(best_schedule.storage_power, 1);
solar_power = grid_data.renewable.solar_capacity * grid_data.renewable.solar_profile;
wind_power = grid_data.renewable.wind_capacity * grid_data.renewable.wind_profile;
total_renewable = solar_power + wind_power;

plot(1:24, total_gen, 'b-', 'LineWidth', 2);
hold on;
plot(1:24, total_load, 'r-', 'LineWidth', 2);
plot(1:24, total_renewable, 'g-', 'LineWidth', 2);
plot(1:24, total_gen + total_renewable - total_storage, 'k--', 'LineWidth', 2);
title('功率平衡分析', 'FontSize', 12, 'FontWeight', 'bold');
xlabel('时间 (小时)', 'FontSize', 10);
ylabel('功率 (MW)', 'FontSize', 10);
grid on;
legend('总发电量', '总负荷', '新能源出力', '净发电量', 'Location', 'best');

% 子图6：成本分解饼图
subplot(2, 3, 6);
cost_labels = {'发电成本', '网络损耗', '碳排放', '储能成本', '约束惩罚'};
cost_values = [cost_breakdown.generation_cost, cost_breakdown.network_loss, ...
               cost_breakdown.carbon_cost, cost_breakdown.storage_cost, ...
               cost_breakdown.violation_penalty];
cost_percentages = 100 * cost_values / final_cost;

pie(cost_percentages, cost_labels);
title('成本分解', 'FontSize', 12, 'FontWeight', 'bold');

% 调整子图间距
sgtitle('智能电网调度优化结果', 'FontSize', 16, 'FontWeight', 'bold');

%% 6. 性能统计
fprintf('\n6. 性能统计：\n');
fprintf('- 总运行时间: %.2f 秒\n', total_time);
fprintf('- 平均每次评估时间: %.4f 秒\n', total_time/algorithm_params.max_evaluations);
fprintf('- 收敛速度: 在%d次评估后达到最优解\n', optimization_history(end, 1));
fprintf('- 成本改善: 从初始%.2f元优化到%.2f元\n', ...
    optimization_history(1, 2), optimization_history(end, 2));

%% 7. 保存结果
fprintf('\n7. 保存结果...\n');
save('grid_optimization_results.mat', 'best_schedule', 'optimization_history', ...
     'total_time', 'grid_data', 'algorithm_params', 'cost_breakdown');

% 生成结果报告
fid = fopen('grid_optimization_report.txt', 'w');
fprintf(fid, '智能电网调度优化报告\n');
fprintf(fid, '====================\n\n');
fprintf(fid, '算法参数：\n');
fprintf(fid, '- 种群大小: %d\n', algorithm_params.population_size);
fprintf(fid, '- 初始样本数: %d\n', algorithm_params.initial_samples);
fprintf(fid, '- 子问题数量: %d\n', algorithm_params.sub_problems);
fprintf(fid, '- 最大评估次数: %d\n', algorithm_params.max_evaluations);
fprintf(fid, '\n优化结果：\n');
fprintf(fid, '- 总运行时间: %.2f 秒\n', total_time);
fprintf(fid, '- 最终成本: %.2f 元\n', final_cost);
fprintf(fid, '- 发电成本: %.2f 元\n', cost_breakdown.generation_cost);
fprintf(fid, '- 网络损耗: %.2f 元\n', cost_breakdown.network_loss);
fprintf(fid, '- 碳排放成本: %.2f 元\n', cost_breakdown.carbon_cost);
fprintf(fid, '- 储能成本: %.2f 元\n', cost_breakdown.storage_cost);
fprintf(fid, '- 约束违反惩罚: %.2f 元\n', cost_breakdown.violation_penalty);
fclose(fid);

fprintf('测试完成！结果已保存到以下文件：\n');
fprintf('- grid_optimization_results.mat (MATLAB数据文件)\n');
fprintf('- grid_optimization_report.txt (文本报告)\n');

% 恢复警告设置
warning('on', 'MATLAB:nearlySingularMatrix');
warning('on', 'MATLAB:singularMatrix');

